Ellipsoid fitting#
Napari-stress implements several algorithms for ellipse fitting. They all have in common that they return a napari vectors layer which represents the three major axes of the fitted ellipsoid. Napari-stress provides some further functionality to project the input pointcloud onto the surface of the fitted ellipsoid.
Napari-stress implementation Least squares ellipsoid fitting
Vedo implementation Statistical ellipsoid fitting
Quantification Quantify fit quality
Mean curvature Measure mean curvature on ellipsoid
import napari_stress
from napari_stress import approximation, measurements, vectors
import napari
import numpy as np
import matplotlib.pyplot as plt
napari.manifest -> 'napari' could not be imported: Could not find file 'builtins.yaml' in module 'napari'
viewer = napari.Viewer(ndisplay=3)
napari.manifest -> 'napari' could not be imported: Could not find file 'builtins.yaml' in module 'napari'
Napari-stress implementation #
This is a least-squares approach at ellipse fitting.
pointcloud = napari_stress.get_droplet_point_cloud()[0][0][:, 1:]
expander = approximation.EllipsoidExpander(get_measurements=True)
expander.fit(pointcloud)
ellipsoid_stress = expander.coefficients_
viewer.layers.clear()
viewer.add_points(pointcloud, size=0.5, face_color='orange')
viewer.add_vectors(ellipsoid_stress, edge_width=1, edge_color='magenta')
napari.utils.nbscreenshot(viewer)
To display, where the initial input points would fall on the surface of the fitted ellipse, use the expand_points_on_ellipse() function:
fitted_points_stress = expander.fit_expand(pointcloud)
viewer.add_points(fitted_points_stress, size=0.5, face_color='cyan')
napari.utils.nbscreenshot(viewer)
Vedo implementation #
This function re-implements the respective function from the vedo library. It applies a PCA-algorithm to a pointcloud to retrieve the major and minor axises of an ellipsoid, that comprises a set fraction of points within its volumne. The inside_fraction parameter controls how many points of the pointcloud will be located within the volume of the determined ellipsoid.
ellipsoid_vedo = napari_stress.fit_ellipsoid_to_pointcloud_vectors(pointcloud, inside_fraction=0.675)
expander.coefficients_ = ellipsoid_vedo
fitted_points_vedo = expander.expand(pointcloud)
viewer.layers.clear()
viewer.add_points(pointcloud, size=0.5, face_color='orange')
viewer.add_vectors(ellipsoid_vedo, edge_width=1, edge_color='magenta')
viewer.add_points(fitted_points_vedo, size=0.5, face_color='cyan')
napari.utils.nbscreenshot(viewer)
Fit quality quantification #
Lastly, if you wanted to quantify the fit remainder (i.e., the distance between the input and the fitted points), you can do this directly using the provided EllipsoidExpander object’:
expander.properties
| euclidian_distance | |
|---|---|
| 0 | 0.873050 |
| 1 | 0.936172 |
| 2 | 0.923862 |
| 3 | 0.905663 |
| 4 | 0.832209 |
| ... | ... |
| 1019 | 0.754312 |
| 1020 | 0.805818 |
| 1021 | 1.071268 |
| 1022 | 0.776711 |
| 1023 | 0.915282 |
1024 rows × 1 columns
You can calculate the length of these vectors using numpy:
distance_vectors = vectors.pairwise_point_distances(pointcloud, fitted_points_vedo)
viewer.layers.clear()
viewer.add_vectors(distance_vectors, edge_width=0.3, edge_color='euclidian_distance', features=expander.properties)
napari.utils.nbscreenshot(viewer)
Mean curvature on ellipsoid #
Lastly, the curvature on the surface of an ellipsoid can be calculated esily with measurements.curvature_on_ellipsoid. Note: The bject returned from this function is of type LayerDataTuple, the structure of which is explained in the docstring. Use help(measurements.curvature_on_ellipsoid) to show.
result = measurements.curvature_on_ellipsoid(ellipsoid_stress, fitted_points_stress)
You can print the global curvatures on the ellipsoid: \(H_e\) (average mean curvature \(H_{e, 0}\) and maximal/minimal mean curvatures \(H_{e, M}\) and \(H_{e, N}\))
for key in result[1]['metadata'].keys():
print(key, ': ', result[1]['metadata'][key])
H0_ellipsoid : 0.0663293556360117
H_ellipsoid_major_medial_minor : [0.08211624 0.06563347 0.05525324]
frame : 0
Plot a histogram of curvatures…
fig, ax = plt.subplots()
hist = ax.hist(result[1]['features']['mean_curvature'], 50)
ax.set_ylabel('Occurrences [#]')
ax.set_xlabel('Mean curvature')
Text(0.5, 0, 'Mean curvature')
…or visualize curvature in the napari viewer:
viewer.layers.clear()
viewer.add_points(result[0], **result[1])
napari.utils.nbscreenshot(viewer, canvas_only=True)
Normals on ellispoid#
You can also calculate the normals on the ellipsoid:
normals = approximation.normals_on_ellipsoid(fitted_points_stress)
viewer.layers.clear()
viewer.add_vectors(normals, edge_width=0.3)
napari.utils.nbscreenshot(viewer)
